2.3. Ellipsoid fitting#

Napari-stress implements several algorithms for ellipse fitting. They all have in common that they return a napari vectors layer which represents the three major axes of the fitted ellipsoid. Napari-stress provides some further functionality to project the input pointcloud onto the surface of the fitted ellipsoid.

import napari_stress
from napari_stress import approximation, measurements
import napari
import numpy as np

import matplotlib.pyplot as plt

---------------------------------------------------------------------------
ModuleNotFoundError                       Traceback (most recent call last)
/tmp/ipykernel_1776/1724537087.py in <module>
----> 1 import napari_stress
      2 from napari_stress import approximation, measurements
      3 import napari
      4 import numpy as np
      5 

ModuleNotFoundError: No module named 'napari_stress'

2.3.1. Napari-stress implementation #

This is a least-squares approach at ellipse fitting.

pointcloud = napari_stress.get_droplet_point_cloud()[0][0][:, 1:]
ellipsoid_stress = approximation.least_squares_ellipsoid(pointcloud)

viewer = napari.Viewer(ndisplay=3)
viewer.add_points(pointcloud, size=0.5, face_color='orange')
viewer.add_vectors(ellipsoid_stress, edge_width=1, edge_color='magenta')
napari.utils.nbscreenshot(viewer)

WARNING: DirectWrite: CreateFontFaceFromHDC() failed (Indicates an error in an input file such as a font file.) for QFontDef(Family="8514oem", pointsize=12, pixelsize=20, styleHint=5, weight=50, stretch=100, hintingPreference=0) LOGFONT("8514oem", lfWidth=0, lfHeight=-20) dpi=192

Assistant skips harvesting pyclesperanto as it's not installed.

To display, where the initial input points would fall on the surface of the fitted ellipse, use the expand_points_on_ellipse() function:

fitted_points_stress = approximation.expand_points_on_ellipse(ellipsoid_stress, pointcloud)
viewer.add_points(fitted_points_stress, size=0.5, face_color='cyan')
napari.utils.nbscreenshot(viewer)

2.3.2. Vedo implementation #

This function re-implements the respective function from the vedo library. It applies a PCA-algorithm to a pointcloud to retrieve the major and minor axises of an ellipsoid, that comprises a set fraction of points within its volumne. The inside_fraction parameter controls how many points of the pointcloud will be located within the volume of the determined ellipsoid.

viewer2 = napari.Viewer(ndisplay=3)
viewer2.add_points(pointcloud, size=0.5, face_color='orange')

<Points layer 'pointcloud' at 0x1cd0e3a8af0>

ellipsoid_vedo = napari_stress.fit_ellipsoid_to_pointcloud_vectors(pointcloud, inside_fraction=0.675)
fitted_points_vedo = approximation.expand_points_on_ellipse(ellipsoid_vedo, pointcloud)

viewer2.add_vectors(ellipsoid_vedo, edge_width=1, edge_color='magenta')
viewer2.add_points(fitted_points_vedo, size=0.5, face_color='cyan')
napari.utils.nbscreenshot(viewer2)

2.3.3. Fit quality quantification #

Lastly, if you wanted to quantify the fit remainder (i.e., the distance between the input and the fitted points), you can do this with approximation.pairwise_point_distances():

residues = approximation.pairwise_point_distances(pointcloud, fitted_points_stress)

You can calculate the length of these vectors using numpy:

lengths = np.linalg.norm(residues[:, 1], axis=1)

viewer.add_vectors(residues, edge_width=1, features = {'lengths': lengths}, edge_color = 'lengths', edge_colormap = 'twilight')
napari.utils.nbscreenshot(viewer)

2.3.4. Mean curvature on ellipsoid #

Lastly, the curvature on the surface of an ellipsoid can be calculated esily with measurements.curvature_on_ellipsoid:

data, features, metadata = measurements.curvature_on_ellipsoid(ellipsoid_stress, fitted_points_stress)

Plot a histogram of curvatures…

fig, ax = plt.subplots()
hist = ax.hist(features['mean_curvature'], 50)
ax.set_ylabel('Occurrences [#]')
ax.set_xlabel('Mean curvature')

Text(0.5, 0, 'Mean curvature')
../../_images/demo_fit_ellipsoid_19_1.png

…or visualize curvature in the napari viewer:

# make other layers invisble
for layer in viewer.layers:
    layer.visible = False

viewer.add_points(data, size=0.5, features=features, face_color='mean_curvature')
napari.utils.nbscreenshot(viewer, canvas_only=True)